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In this paper we introduce a three-dimensional version of the Mercedes-Benz model to describe 
water molecules. In this model van der Waals interactions and hydrogen bonds are given explicitly 
through a Lennard- Jones potential and a Gaussian orientation-dependent terms, respectively. At low 
temperature the model freezes forming Ice-I and it reproduces the main peaks of the experimental 
radial distribution function of water. In addition to these structural properties, the model also 
captures the thermodynamical anomalies of water: the anomalous density profile, the negative 
thermal expansivity, the large heat capacity and the minimum in the isothermal compressibility. 



I. INTRODUCTION 

Water is the most important fluid on earth. It covers 
two thirds of the planet's surface and controls its cli- 
mate. Most importantly, water is necessary for carbon- 
based organic life being the solvent in most in vivo chem- 
ical reactions. Its unique hydration properties drive bi- 
ological macromolecules towards their three-dimensional 
structure, thus accounting for their function in living or- 
ganisms [T]. Water exhibits anomalous properties that 
affect life at a larger scale. For example, mammals bene- 
flt from the large latent heat of water to cool them down 
through sweating, while water's large heat capacity pre- 
vents local temperature fluctuations, facilitating thermal 
regulation of organisms. 

These anomalous properties result from a competition 
between isotropic van der Waals interactions and highly 
directional hydrogen bonding (H-bond). A large number 
of models of varying complexity have been developed and 
analyzed to model water's extraordinary properties, for 
reviews see e.g. Refs. [3J [31 HI [S] , but none of the current 
models can correctly reproduce all physical properties of 
water. Those model are typically calibrated against ex- 
perimental data, for example the radial distribution func- 
tion (RDF) at ambient conditions OE], or the temper- 
ature of maximum density [8], i.e. T = 3.98°C. While 
there is no guarantee that a model optimized to repro- 
duce a given property is able to account for others, adding 
details increases its quantitative accuracy. For example, 
TIP5P, which describes water through 5 interacting sites, 
is typically more accurate [5] than models with 3 or 4 
interacting sites. The addition of each interacting site, 
however, makes the model considerably much more de- 
manding computationally. This is an undesirable feature 
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since a large number of water molecules is required to hy- 
drate even the smallest peptides resulting in a high com- 
putational cost. Thus, simple models, such as SPC [S] 
and TIP3P [TOl E] , are the most used ones in computa- 
tional studies of biologically motivated systems. In ad- 
dition, new models and improvements appear frequently 
in literature, see e.g., Refs. [T^l [I3J HH [TS] and references 
in them. 

Coarse-grained models have also been developed and 
used to study the emergence of water's anomalous prop- 
erties from its atomic constituents. Both lattice [T6j [171 
[T5] and continuous models [ini [Ml IHl HI] have been ap- 
plied. Current coarse-grained models cannot, however, 
be easily used to study hydration of macromolecules, be- 
cause they can not reproduce the structure of liquid wa- 
ter which is essential in studies of biological systems and 
molecules [T]. A proper structural description is required 
since hydration and, in particular, the hydrophobic effect, 
which is the main driving force for protein folding [23. 24] , 
depend on the amount of structural order close to the hy- 
drated molecule versus the amount of order in bulk wa- 
ter. A simplified model that would account for both ther- 
modynamical and structural properties of water, would 
therefore be highly beneficial in studies related to the hy- 
drophobic effect, protein folding and macromolecules in 
general. 

The main purpose of this work is to introduce a sim- 
ple but realistic model that reproduces both the main 
structural and thermodynamic properties of water. To 
this end, we extend the two dimensional (2D) Mercedes- 
Benz (MB) model [25] to 3D. In 2D, the MB model has 
already provided insights into several properties of wa- 
ter: its anomalous thermodynamical behavior |19j , hy- 
dration of non-polar solutes [26], ion solvation [27 , cold 
denaturation of proteins [28', and the properties of dif- 
ferent amino acids [29 . Despite this success, there are 
several mechanisms which cannot be studied in two di- 
mensions and an extension to 3D is needed. 

We show that a previously proposed framework for the 
3D MB model [30] does not reproduce the thermodynam- 
ical anomalies of water. Here, we extend the model to 
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FIG. 1: (Color online) Schematic representation of two 
Mercedes-Benz molecules and important vectors defining 
their interaction. 



overcome that problem by making H-bonding dependent 
on the local environment of atoms by penalizing compact 
configurations in favor of open-packed ones. With this 
implementation, structural and thermodynamical prop- 
erties of are recovered qualitatively. We would like to 
emphasize that our goal is not to replace atomistically 
accurate models, such as TIP5P, but rather to provide 
an alternative for qualitative studies involving water. 

The rest of this paper is organized as follows: next, 
we introduce and discuss the existing 3D MB model, and 
propose a correction that makes the model suitable to 
describe thermodynamical and structural properties of 
water. In the same section, we present the Monte-Carlo 
scheme and the cooling protocol used in this work. We 
present our results and comparison to experimental re- 



sults in Sec. Ill In the section entitled results, experimen- 
tal data are compared qualitatively to our simulations. 
Finally, we present our conclusions and a discussion in 
Sec.lTTl 



II. THE MODEL 

A. Mercedes Benz model 

In the 3D MB model, water molecules interact explic- 
itly through two types of empirical potentials: H-bonds 
and van der Waals. H-bonds are directional and account 
for the tetrahedral structure of water which is described 
by four arms separated from each other by angles of 
109.47°, see Fig. [T] The energy of H-bonds is minimized 
whenever arms of adjacent molecules point towards each 
other. Mathematically if Xi represents the position of 
the i"^ particle and its four unitary arms, which are de- 
noted by ifc (with A: = 1, 2, 3, 4), then H-bond interaction 
between molecules i and j can be written as 



k,l = l 



(1) 



where 

G{ikfij - l,ae)G{jifij - l,crg) (2) 

and G{x, a) is an unnormalized Gaussian function 

G{x,a) ^exp[-x^/2a% (3) 

The above mathematical description ensures that the in- 
tensity of a H-bond is maximized whenever the arms 
of neighboring molecules are aligned with the vector r^j 
joining their centers of mass and whenever their distance 
is equal to Rhb- 

The spherically symmetric van der Waals interactions 
are approximated by a Lennard-Jones potential: 



ULj{nj) = 4e 



LJ 



12 



(4) 



where clj describes the strength of the interaction and 
CT^j is the particle diameter. Then, the total energy de- 
scribing two MB particles is given by 



U{X,,X,) - UlM,) + Uhb{X,,X,). 



(5) 



Bizjak et al. [30J studied this model using the following 
set of parameters: euB = —1, = l/^Se/fs, Rhb = 1, 
CT^j — 0.7, uji — a0 — 0.085. They assumed a diamond 
structure for the model's ground state which is the config- 
uration taken by oxygen atoms when water forms cubic- 
ice, i.e. Ice-Ic. When tested against simulation, how- 
ever, this assumption fails and the model can be shown 
to minimize its energy in an Ice- VII configuration: two 
interpenetrating diamond lattices with no H-bonds con- 
necting these lattices. 

Ice- VII appears to be the optimized ground state for 
systems trying to maximize their density within a tetra- 
hedral symmetry. It is therefore natural that the 3D 
MB model of Bizjak et al. |30], whose H-bond term im- 
poses a tetrahedral configuration and the van der Waals 
term favors compact conformations, has this structure 
as its ground state. However Ice- VII is not the desired 
ground state for models of water at ambient pressure such 
that 3D lattice models for this material have an explicit 
energetic term penalizing compact configurations of this 
type [18, 31, 32]. 

Figure [2][a) shows a typical Ice-VII ground state ob- 
tained by quenching a system of 256 MB particles inter- 
acting through the framework of Bizjak et al.. Details 
about the simulation method and the cooling procedure 
will be described later in this section. Figure [2]jb) shows 
the density obtained along quenching. The liquid phase 
has a higher density than ice ~ as in real water. How- 
ever, the model does not reproduce the density anomaly 
of water: the liquid phase does not show a temperature 
of maximum density below which the density decreases. 
As a result, the thermal expansion coefficient is never 
found to be negative in the liquid phase (Fig. |2jb)). In 
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FIG. 2: (Color online) (a) Ground state of the Mercedes- 
Benz model: Ice- VII. As a guide to the eye two hexagons, 
representing the two interpenetrating diamond structure of 
ice- VII, are drawn, (b - left axis) Dependence of the density 
(g/cm^) on temperature (in units of e£rs). (b - right axis) 
Dependence of the coefficient of thermal expansion (e^^) on 
temperature. A constant pressure of 0.2 in units of ens/RjiB 
was used. 



addition, the model does not reproduce the structure of 
liquid water (see Fig. 3 in Ref. [30 ): the simulated RDF 
has a non-realistic peak at a distance corresponding to 
the van der Waals radius. 

Despite the problems cited above, the 3D MB model 
remains an attractive coarse-grained model for water. It 
does not require calculation of charges, which enables 
longer simulation times desperately needed in studies 
of macromolecules. It also holds the promise of being 
able to provide a qualitatively accurate description of the 
structure of water due to its tetrahedral nature |33] and 
water's thermodynamical properties since the model ex- 
hibits both open and close packed structures required to 
described water's anomalous behavior PU]- Next, we de- 



scribe how the model of Bizjak et al. [30] can be improved 
to better describe the structure and thermodynamics of 
water. 



B. Corrections to the Mercedes-Benz model 

To resolve the above problems, we introduce a term 
that depends on the local environment of particles. Our 
approach is inspired by TersofF-like potentials for covalent 
materials 34J. This term penalizes H-bonds which are 
formed in crowded environments through the factor 

r 1, z/z, < 4 

where Zi is the coordination of atom i, computed as 
■^i ~ '^k^i firik) with the cut-off function defined 
by W:- 

{1, r<R-D 
isin(f(r-i?)/i:>), R- D <r < R + D 
0, r>R + D 

(7) 

where R and D are chosen as to include the first- neighbor 
shell only. Note that /(r) decreases continuously from 1 
to in the range R — D<r<R+D. 

The energy of H-bonds corrected through Eq. ^ be- 
comes: 

4 

u^HBiXr,x,) = biz,) J2 ul}sinj,il,fi). (8) 

k,l = l 

This equation penalizes H-bonds whenever interacting 
molecules have more than four neighbors. This inhibits 
the formation of compact tetrahedral phases, e.g. Ice- 
VII, and favors open-packed tetrahedral phases such as 
Ice-I. 

In order to ensure that H-bonds favor chair-like con- 
figurations required for diamond structure, we also add 
a standard potential with three-fold symmetry for dihe- 
dral angles. This potential adds an energetic cost to 
the H-bond between arms m of molecule i and arm n 
of molecules j, if the dihedral angle formed by the other 
arms of these molecules is not 60°: 

C/^'"(X„ X,-) = |c/5'S(^.j, X (9) 

^(l+cos(30feO), 

1^71 

where is the strength of the interaction. The term 
U^g{rij,i„i, jn)b{zi) ensures that the penalty is propor- 
tional to the strength of the H-bond. The dihedral angle 
(j)ki describes how the arm k of molecule i aligns with the 
arm I of molecule j along the vector joining the center 
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of mass of these two molecules. Thus, the total dihedral 
energy between molecules i and j is 



(10) 



Note that because of the dependence on the lo- 
cal environment, Ufjg{Xi,Xj) Ufjg{Xj, Xi) and 
C/^™(X,,Xj) 7^ [/J"(X,-,X,). This asymmetry has no 
physical implications since Ufj^ and possess all the 
invariance properties required for a potential [34_ . 

We can now write the total potential energy between 
two water molecules as 

E{X,,Xj) = UlMj) + U'hb{X,,X,) + U4X,,Xj). 

(11) 

This model has 10 parameters which were chosen such as 
to account for a semi-quantitative agreement of the den- 
sity profile with experiment. We proceeded in two steps 
to adjust these parameters. First, we chose the values for 
these parameters such as to produce a density in g/cm'^ 
[37j that is comparable to experimental values, i.e. about 
1 g/cm"^ for the liquid phase and 0.93 g/cm"^ for the ice 
phase. Only under this condition can the structure of 
the model be qualitatively similar to real water. Then, 
we adjusted the parameters such as to obtain a density 
that is a concave function of temperature with its max- 
imum close to the freezing point. This second condition 
is the minimal requirement for describing the anomalous 
properties of water. 

The set of parameters calibrated according to the 
above procedure is given here in reduced units. We re- 
port energies and distances in terms of the binding en- 
ergy lensl and equilibrium distance Rhb of the H-bond. 
In these units, the three binding energies describing the 
system are ens = —1, clj = 0.05 and = 0.01. The 
two distances are Rhb = 1 and aLj — 1.04/2^/^. The 
two terms controlling H-bond interaction are ~ 0.1 
and (79 — 0.08, and the three parameters controlling the 
penalty of crowded environments are v = 0.5, i? = 1.3 
and D = 0.2. In this work, temperature is given in units 
of leHsl/fcs, where Boltzmann's constant fcs is set to 
unity. Pressure is given in units of \tHB\l R?hb- 

While adjusting the parameters, we found that the be- 
havior of the system is robust upon changing the vari- 
ables characterizing crowded environments. It is, how- 
ever, sensitive to the ratio between the binding energy of 
the van der Waals interaction and the binding energy of 
the H-bond. This ratio controls the interplay of forces 
leading to an environment where MB molecules are ra- 
dially surrounded by their first- neighbors, and forces fa- 
voring a tetrahedral distribution of the first-neighbors. 
The latter favors a high density configuration while the 
former accounts for a low density one. As opposed to the 
2D MB model, we kept the equilibrium distance of the 
van der Waals interaction comparable to the equilibrium 
distance of the H-bond such as to avoid artificial peaks 
in the RDF (30j . 



C. Numerical simulation method 

For numerical simulations, we use the isothermal- 
isobaric (NPT) ensemble to study the thermodynamical 
properties of a system made of TV = 256 MB particles. 
A Monte-Carlo scheme is used where, at each step, an 
attempt is made to displace the center of mass and the 
orientation of particles randomly by a quantity ^R^^x 
and 0.125 rad, respectively. The maximum translational 
displacement is chosen such as to give an acceptance ratio 
of 50%. Periodic boundary conditions are used to mimic 
an infinite system and at every 5 Monte Carlo sweeps, an 
attempt to rescale the size of the box is made (1 Monte 
Carlo sweep is equivalent to X attempted steps). 

To obtain thermodynamical data throughout the de- 
sired range of temperatures, the initial configuration of 
the system is chosen randomly and equilibrated at the 
highest temperature (T = 0.17) for 5 x 10'* sweeps, af- 
ter which statistics are gathered for the same amount of 
time. Then, the system is cooled down by AT — 0.002 
and a similar cycle of equilibration/data gathering is per- 
formed. This cooling procedure is repeated until the 
lowest temperature, i.e. T = 0.11 is reached. At the 
transition temperature an additional cycle of equilibra- 
tion/statistics gathering ensured that the system was 
equilibrated properly. For all the pressures studied here, 
this protocol was repeated for 10 samples difi^ering by 
the initial condition. All the quantities reported are the 
average over those 10 samples and, whenever relevant, 
the root-mean-square of this average is also shown as the 
error-bar. 

The quantities computed during the simulations were 
the average potential energy per particle _B, the volume 
per particle V , the heat capacity Cp, the compressibility 
and the thermal expansion coefficient ap. The last 
three quantities are computed mathematically from the 
standard fluctuation relations: 



(12) 



T{V) ' 
, ^ {VH)-{V){H) 

where H corresponds to the enthalpy of the system. As 
for the other quantities computed during the simulation, 
these response functions will be given in reduced units. 
Thus, Cp will be reported in dimensionless units, in 
terms of R^^^/tHB and a*p in units of e^^. 



III. RESULTS 

In Fig. [3] we provide a qualitative comparison between 
the properties of bulk water (left panels) and the MB 
model at P = 0.2 (right panels). The behavior of the 
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FIG. 3: Thermodynamical properties of water. Left: Experiments - data obtained from references [35] . Units are: p (g/cm^), 
Op (K~^), Cp (J / g / K) and k.t (bar~^). Right: Simulations. Units are: Op k.t {R'hb /^hb) and Cp (dimensionless 

units). 
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MB model follows the trends of water quite accurately: 
the anomalous density profile (panels on the first row), 
the negative thermal expansivity (second row), the min- 
imum in the isothermal compressibility (third row) and 
the large heat capacity (fourth row). 

At ambient pressure, water freezes into an open packed 
configuration called hexagonal-ice, i.e. Ice-lh- This struc- 
ture is held together by H-bonds which break when ice 
melts. At this transition, water molecules fill part of 
the empty spaces, assuming a higher density. In the liq- 
uid phase close to the melting temperature, a few open- 
packed configurations persist - held together by H-bonds. 
As the system is heated up, those bonds melt gradu- 
ally removing empty spaces and increasing the density of 
the system. This reduction of empty spaces occurs until 
the temperature of maximum density is reached. At this 
point thermal fiuctuations decrease the density of the sys- 
tem with increasing temperature. This behavior has been 
measured experimentally (Fig. |3ja)) and is captured by 
the MB model (Fig. [sj^b)): abrupt increase of the density 
at the melting transition and concave temperature de- 
pendence for the density of water with a maximum close 
to the melting transition. 



The thermal expansion coefficient is proportional to 
the derivative of the volume with respect to tempera- 
ture ap = l/V{dV/dT)p. As for most materials, ap 
decreases upon coohng (Fig. |3jc)) - indicating that the 
volume of water decreases with temperature. It becomes 
zero at the temperature of maximum density and neg- 
ative close to the freezing point. This unusual negative 
expansivity is reproduced in the model (Fig. |3jd)) and 
reflects the unusual behavior of water to expand upon 
cooling below the temperature of maximum density. 

The isothermal compressibility measures the tendency 
of a system to change its volume when the applied pres- 
sure is varied: kt = —l/V{dV/dP)T- For a typical ma- 
terial, kt decreases upon cooling since it is related to 
density fluctuations whose amplitude becomes smaller as 
temperature decreases. This is in contrast with the be- 
havior of water (Fig.jsje)). The compressibility of water 
is a convex function of temperature and has a minimum. 
This anomalous behavior can be explained by noticing 
that the compressibility is lower for highly packed system 
than for loosely packed ones since highly packed systems 
are less susceptible to rearrange their conformation when 
subjected to a pressure change. Thus, kt correlates with 
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FIG. 4: Schematical representation of the phase diagram of a 
simple one-component substance (a) and water (b). Arrows 
indicate that pressure freezes a typical liquid but pressure 
melts ice. (c) Simulated density of the model for different 
values of pressure. Arrow indicates the shift of the freezing 
temperature to lower values as pressure increases. 



the volume of the system [T^]. Now, since the volume of 
water is a convex function of temperature, kt is also con- 
vex with respect to temperature for water ~ see Fig.[3][e). 
Fig. [sf^f ) shows that the simulated compressibiHty is also 
a convex function of temperature, although the curva- 
ture is not very pronounced and its minimum is not as 
pronounced as in experiments. 

Heat capacity, which measures the capacity of a sys- 
tem to store thermal energy (Cp = {dH/dT)p), is much 
higher in water than in ice - see Fig.[3|h). This has been 
explained by the multiple energy storage mechanisms of 
water as the breakage of van der Waals interactions and 
H-bonds. The heat capacity of the model presents a much 
higher variability than real water: close to the transition, 
Cp is much higher than ice and this quantity decreases 
fast, reaching the same value as ice at about T = 0.15. 

In Fig.|4|^a), we illustrate schematically the coexistence 
lines of the solid, liquid and vapor phases of a simple 
material. At any point along those lines, the free energies 
of the adjacent phases are equivalent and the Clausius- 
Clapeyron equation is obtained by equating them: 




\dTJ coex TAv 



Ah 



(13) 



Since Aft, < and Au < for the liquid to solid transi- 
tion of typical materials, {dP/dT)coex is positive. As a 
result of this positive slope, a typical liquid freezes when 
pressure is applied to it - as illustrated by the arrow on 
Fig. Qa). On the other side, since water expands upon 
freezing, Av > while the enthalpy difference remains 
negative (ice has a lower enthalpy compared to water). 
Thus the coexistence line of the liquid-solid transition has 
a negative slope, i.e. {dP/dT)coex < 0. This is illustrated 
in Fig. Qb) and leads to the melting of ice when pres- 
sure is applied to it. In Fig.Qc) we show that the model 
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FIG. 5: (Color online) (a) Ground state of the model: Ice-I 
where oxygen atoms occupy positions on a diamond-like lat- 
tice, (b) Radial distribution function of the model at different 
temperatures compared to experimental data of water at 298 
K. 



reproduces this anomalous behavior of water. The simu- 
lated dependence of the density on temperature is shown 
for different values of pressure. The freezing tempera- 
ture shifts to lower values as pressure increases, implying 

idP/dT),oex < 0. 

Figure |5ja) shows the structure obtained by freez- 
ing the MB model. In this configuration, the cen- 
ter of a MB molecule occupies the sites of a diamond- 
like structure and its arms point to its four neighbors. 
Note that without a penalty term (Eq. ^) the empty 
spaces found in Ice-I can be the stage for the forma- 
tion of another tetrahedral-like lattice. Thus, this term 
efficiently shifts the energy of those compact configu- 
rations and, in particular Ice- VII, in favor of Ice-I. In 
Fig. [sjjb), the experimental RDF [7^ is compared to the 
ones of the model at different temperature. The second 
peak of the RDF, commonly referred to as tetrahedral 
peak |36j , is a fingerprint of the tetrahedral geometry of 
water. It occurs at a distance given by the cosine rule, 

= 2RjjB-2RHB cos(109.4°) « 1.6, much smaller than 
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for a simple liquid [351. When compared to experiment, 
the model's RDF has slightly less structure but it peaks 
at the same position as the experiment - indicating that 
the average structure of the model agrees well with the 
experiment. 

IV. CONCLUSION 

In this work, we have constructed a simple but realistic 
model for water based on the Mercedez-Benz approach 
P^ . At low temperature the model freezes forming Ice- 
I and it reproduces the main peaks of the experimen- 
tal RDF of bulk water. In addition to these structural 
properties, the model reproduces the density anomaly of 
water: ice has a lower density than water and the den- 
sity of water is a concave function of temperature, with 
a maximum close to the freezing point. Also, the slope 
of the solid-liquid coexistence curve is also found to be 
negative, in agreement with experiments. 

In the 2D MB model, the H-bond interaction favors 
environments having 3 first-neighbors at a distance Rh b 
which competes with the van der Waals interaction that 
favors 6 neighbors at a distance 0.7 Rhb- This competi- 
tion is the underlying physics of the model that accounts 
for the density anomaly of bulk water. In the 3D MB 
model, H-bond and van der Waals interactions have the 
same equilibrium distance and the density anomaly re- 



sults from an energy penalty for crowded environments. 
Without this penalty, the system would solidify into a 
compact Ice- VII configuration. With the penalty term, 
Ice- VII conformations compete with an open packed 
diamond-like structure. The competition between these 
interactions is the underlying mechanism that leads to 
the density anomaly of the system. 

The MB model for water is based on local interactions 
which are much faster to compute than usual models that 
uses long-range Coulomb forces. We believe that this 
model will provide new insights into water mechanisms 
related to molecular hydration. In particular, investiga- 
tions of the hydrophobic effect are being undertaken with 
this model. 



Acknowledgements 

C.L.D. would like to thank Razvan Nistor and Marco 
Aurelio Alves Barbosa for insightful discussions. He 
would also like to thank Alvarro Ferraz Filho and Sil- 
vio Quczado for kindly hosting his stay at the Interna- 
tional Centre of Condensed Matter Physics (ICCMP) in 
Brasilia, Brazil, where part of this work was completed. 
We would like to thank SharcNet (www.sharcnet.ca) 
for computing resources. M.K. has been supported by 
NSERC of Canada and T.A-N. by the Academy of Fin- 
land through its COMP CoE and TransPoly grants. 



[1] Martin Chaplin. Do we underestimate the importance 
of water in cell biology? Nature Reviews Molecular Cell 
Biology, 7(ll):861-866, September 2006. 

[2] William L. Jorgensen, Jayaraman Chandrasekhar, Jef- 
fry D. Madura, Roger W. Impey, and Michael L. Klein. 
Comparison of simple potential functions for simulat- 
ing liquid water. The Journal of Chemical Physics, 
79(2):926-935, 1983. 

[3] I. Nezbeda. Simple short-ranged models of water and 
their application, a review. Journal of Molecular Ltqutds, 
73-74:317-336, November 1997. 

[4] Bertrand Guillot. A reappraisal of what we have learnt 
during three decades of computer simulations on water. 
Journal of Molecular Liquids, 101(l-3):219-260, Novem- 
ber 2002. 

[5] M. M. Conde C. Vega, J. L. F. Abascal and J. L. 
Aragones. What ice can teach us about water interac- 
tions: a critical comparison of the performance of differ- 
ent water models. Faraday Discuss., 141:251, 2009. 

[6] Jon M. Sorenson, Greg Hura, Robert M. Glaeser, and 
Teresa Head-Gordon. What can x-ray scattering tell us 
about the radial distribution functions of water? J. 
Chem. Phys., 113:9149, 2000. 

[7] A. K. Soper. The radial distribution functions of water 
and ice from 220 to 673 k and at pressures up to 400 
mpa. Chemical Physics, Volume 258:121-137, 2000. 

[8] Michael W. Mahoney and William L. Jorgensen. A five- 
site model for liquid water and the reproduction of the 
density anomaly by rigid, nonpolarizable potential func- 



tions. J. Chem. Phys., 112:8910, 2000. 
[9] H. J. C. Berendsen, J. P. M. Postma, W. F. van Gun- 
steren, and J. Hermans. Interaction models for water in 
relation to protein hydration. In B. Pullman, editor, In- 
termolecular Forces, pages 331-342. Reidel, Dordrecht, 
1981. 

[10] W. L. Jorgensen, J. Chandrasekhar, J. D. Madura, R. W. 
Impey, and M. L. Klein. Comparison of simple potential 
functions for simulating liquid water. J. Chem. Phys., 
79:926-935, 1983. 

[11] Eyal Neria, Stefan Fischer, and Martin Karplus. Simula- 
tion of activation free energies in molecular systems. J. 
Chem. Phys., 105:1902-1921, 1996. 

[12] Daniel J. Price and Brooks. A modified tip3p water po- 
tential for simulation with ewald summation. The Jour- 
nal of Chemical Physics, 121(20):10096-10103, 2004. 

[13] Hans W. Horn, Wilham C. Swope, Jed W. Pitera, Jef- 
fry D. Madura, Thomas J. Dick, Greg L. Hura, and 
Teresa Head-Gordon. Development of an improved four- 
site water model for biomolecular simulations: Tip4p- 
ew. The Journal of Chemical Physics, 120(20):9665- 
9678, 2004. 

[14] Yujie Wu, Harald L. Tepper, and Gregory A. Voth. 
Flexible simple point-charge water model with improved 
liquid-state properties. The Journal of Chemical Physics, 
124(2), 2006. 

[15] Andras Baranyai and Albert Bartok. Classical interac- 
tion model for the water molecule. The Journal of Chem- 
ical Physics, 126(18), 2007. 



8 



[16] Marco Aurlio Alves Barbosa and Vera Bohomoletz Hen- 
riques. Frustration and anomalous behavior in the bell- 
lavis model of hquid water. Phys. Rev. E 77, 051204 
(2008), 77:51204, 2008. 

[17] C. Buzano, E. De Stefanis, A. Pelizzola, and M. Pretti. 
Two-dimensional lattice-fluid model with waterlike 
anomahes. Phys. Rev. E, 69:61502, 2004. 

[18] M. Pretti and C. Buzano. Thermodynamic anomalies in 
a lattice model of water: Solvation properties. J. Chem. 
Phys., 123:24506, 2005. 

[19] Kevin A. T. Silverstein, A. D. J. Haymet, and Ken A. 
Dill. A simple model of water and the hydrophobic effect. 
J. Am. Chem. Soc., 120 (13):3166-3175, 1998. 

[20] Peter H. Poole, Francesco Sciortino, Tor Grande, H. Eu- 
gene Stanley, and C. Austen Angell. Effect of hydrogen 
bonds on the thermodynamic behavior of liquid water. 
Phys. Rev. Lett, 73:1632-1635, 1994. 

[21] A. P. Lyubartsev and A. Laaksonen. Determination of 
effective pair potentials from ab-initio simulations: Ap- 
plication to liquid water. Chem. Phys. Lett., 325:15-21, 
2000. 

[22] A. P. Lyubartsev, M. Karttunen, I. Vattulainen, and 
A. Laaksonen. On coarse-graining by the Inverse Monte 
Carlo method: Dissipative particle dynamics simulations 
made to a precise tool in soft matter modeling. Soft Ma- 
terials, 1:121-137, 2003. 

[23] W. Kauzmann. Adv. Protein Chem., 14:1, 1959. 

[24] Ken A. Dill. Dominant forces in protein folding. Bio- 
chemistry, 29:7133, 1990. 

[25] A. Ben-Naim. Statistical mechanics of "waterlike" parti- 
cles in two dimensions, i. physical model and application 
of the percus yevick equation. The Journal of Chemical 
Physics, 54:3682, 1971. 

[26] Noel T. Southall and Ken A. Dill. Potential of mean force 
between two hydrophobic solutes in water. Biophysical 
Chemistry, 101-102:295-307, 2002. 

[27] Ken A. Dill, Thomas M. Truskett, Vojko Vlachy, and 
Barbara Hribar-Lee. Modeling water, the hydrophobic 
effect, and ion solvation. Annual Review of Biophysics 
and Biomolecular Structure, 34:173, 2005. 

[28] Cristiano L. Dias, Tapio Ala-Nissila, Mikko Karttunen, 
Ilpo Vattulainen, and Martin Grant. Microscopic mech- 
anism for cold denaturation. Physical Review Letters, 
100:118101, 2008. 

[29] Jean-Paul Becker and Olivier Collet. Mercedes benz 
model of neutral amino-acid side chains. .Journal of 
Molecular Structure: THEOCHEM, 774:23-28, 2006. 

[30] A. Bizjak, T. Urbic, V. Vlachy, and K.A. Dill. The 
three-dimensional " mercedes benz" model of water. Acta 
Chimica Slovenica, 54:532-537, 2007. 

[31] C. J. Roberts and Pablo G. Debenedetti. Polyamorphism 
and density anomalies in network- forming fluids: Zeroth- 
and first-order approximations. J. Chem. Phys., 105:658, 
1996. 

[32] G M Bell. Statistical mechanics of water: lattice model 
with directed bonding. J. Phys. C: Solid State Phys, 
5:889, 1972. 

[33] J. D. Bernal and R. H. Fowler. A theory of water and 
ionic solution, with particular reference to hydrogen and 
hydroxyl ions. J. Chem. Phys., 1:515, 1933. 

[34] J. Tersoff. New empirical approach for the structure and 
energy of covalent systems. Phys. Rev. B, 37:6991-7000, 
1988. 

[35] GS Kell. Density, thermal expansivity, and compressibil- 



ity of liquid water from 0. deg. to 150. deg.: Correlations 
and tables for atmospheric pressure and saturation re- 
viewed and expressed on 1968 temperature scale. Journal 
of Chemical and Engineering Data, 20:97, 1975. 

[36] John L. Finney. The water molecule and its interactions: 
the interaction between theory, modelling, and experi- 
ment. Journal of Molecular Liquids, 90:303-312, 2001. 

[37] The density can be computed in g per cm^ by mapping 
Rhb to its experimental value [7j, i.e. Rhb = 2.78 A, 
and using M = 2.992 • 10"^'^ g for the molecular mass of 
water: p= (^^s") 1.45443 g/cra. 



